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Jellium, a simple model of metals, is a standard testing ground for density functionals, both for 
bulk and for surface properties. Earlier tests show that the Tao-Perdew-Staroverov-Scuseria (TPSS) 
nonempirical meta-generalized gradient approximation (meta-GGA) for the exchange-correlation 
energy yields more accurate surface energies than the local spin density (LSD) approximation for 
spin-unpolarized jellium. In this study, work functions and surface energies of a jellium metal in the 
presence of "internal" and external magnetic fields are calculated with LSD, Perdew-Burke-Ernzerhof 
(PBE) GGA, and TPSS meta-GGA and its predecessor, the nearly nonempirical Perdew-Kurth- 
Zupan-Blaha (PKZB) meta-GGA, using self-consistent LSD orbitals and densities. The results 
show that: (i) For normal bulk densities, the surface correlation energy is the same in TPSS as 
in PBE, as it should be since TPSS strives to represent a self-correlation correction to PBE; (ii) 
Normal surface density profiles can be scaled uniformly to the low-density or strong-interaction 
limit, and TPSS provides an estimate for that limit that is consistent with (but probably more 
accurate than) other estimates; (iii) For both normal and low densities, TPSS provides the same 
description of surface magnetism as PBE, suggesting that these approximations may be generally 
equivalent for magnetism. The energies of jellium spheres with up to 106 electrons are calculated 
using density functionals and compared to those obtained with Diffusion Quantum Monte Carlo 
data, including our estimate for the fixed-node correction. Typically, while PBE energies are too 
low for spheres with more than about two electrons, LSD and TPSS are accurate there. We confirm 
that curvature energies are lower in PBE and TPSS than in LSD. Finally we calculate the linear 
response of bulk jellium using these density functionals, and find that not only LSD but also PBE 
GGA and TPSS meta-GGA yield a linear-response in good agreement with that of the Quantum 
Monte Carlo method, for wavevectors of the perturbing external potential up to twice the Fermi 
wavevector. 

PACS numbers: 71.15.Mb,31.15.Ew,71.45.Gm 



I. INTRODUCTION 

Jellium is a simple model of metals. The surface prop- 
erties of jellium can emulate those of real surfaces. In this 
popular model the lattice of ions is replaced by a uniform 
positive charge density. Lang and Kohni made the first 
self-consistent calculation of the surface properties of this 
model using the local spin density (LSD) approximation 2 
for the exchange-correlation (xc) energy within the Kohn- 
Sham density functional formalism 2 ^^. The results, af- 
ter a correction for lattice effects, agreed surprisingly well 
with experiments. This unexpected success initiated the 
application of the density functional theory to surfaces. 
A restricted variational calculation by Mahan^ gave den- 
sity profiles and surface energies very similar to those of 
Lang and Kohn. 

Perdew and Monnier— performed a series of self- 
consistent calculations of the surface properties of real 
and jellium metals within LSD. Kautz and Schwartz^ ex- 
tended the work of Lang and Kohn to the spin-polarized 
case and calculated self-consistently several surface prop- 



erties of jellium polarized by a magnetic field uniform in- 
side and zero outside the edge of the positive background. 

Over the past few decades, the correct surface 
energy of jellium has been controversial. Recent 
wor k 10 ! 11 ' 12 ! 13 ' 14 ! 15 ! 16 resolves most of the controversy, 
placing this energy close to but probably a little higher 
than that of LSD (which benefits from a strong error 
cancellation between exchange and correlation). For a 
review of this and earlier work, see Ref. [Ifjl 

In density functional theory 2 ^^, everything is treated 
exactly except the exchange-correlation (xc) energy, 
which has to be approximated as a functional of the elec- 
tron density. Development of better density functional 
approximations has been the subject of continuing the- 
oretical efforts. LSD is the simplest approximation and 
has been successful in condensed matter physics. How- 
ever, it tends to overbind molecules. An efficient way 
to solve this problem is to introduce the density gra- 
dient Vn as an additional local ingredient to construct 
a gradient-corrected density functional, the generalized 
gradient approximation (GGA) 17 i 18 . With the advent of 
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GGA, density functional theory has become a popular 
method in quantum chemistry as well. Although GGAs 
have achieved significant improvement over LSD for most 
properties and for diverse systems^, they usually un- 
derestimate the surface exchange-correlation energy of a 
spin-unpolarized jellium, for which LSD is much more 
accurat o 11 ' 20 . 

Further systematic improvement over GGAs may be 
made by imposing additional exact conditions without 
losing those already satisfied by GGA. This can be 
done by employing the kinetic energy densities r-f (r) and 
Tj_(r), and/or the Laplacians of the densities V 2 nf (r) and 
V 2 n^(r), as further additional ingredients. Here the ki- 
netic energy density of electrons with spin a (a =f, j) is 
denned as 

occup 

= 2 E Ml 2 , (!) 

i 

where ipf (r) are the occupied Kohn-Sham orbitals, which 
are implicit functionals of the density n a (r). (Atomic 
units h = ni = e 2 = 1 are used throughout unless other- 
wise explicitly stated.) This family of density functional 
approximations is called meta-GG A 21 i 22 ' 23 . 

While the exact form of the universal functional re- 
mains unknown, many exact constraints on this exact 
functional have been uncovered. Thus the more exact 
constraints a density functional satisfies, the closer it is 
to the exact universal functional. Having this in mind, 
Perdew, Kurth, Zupan, and Blaha (PKZB) 22 constructed 
a meta-GGA from first-principles. This nearly non- 
empirical functional satisfies two important constraints 
which can not be satisfied at the GGA level: It nearly 
recovers the known fourth-order gradient expansion 2 ^ of 
the exchange energy in the slowly- varying limit and it is 
free from self-correlation errors. 

PKZB meta-GGA has impressively correcte d 11 ' 20 ^ 
the too-low surface exchange-correlation energy of GGA 
functionals, and has successfully improved upon LSD and 
GGAs in thermochemistry for molecular atomization en- 
ergies 25 . The defect of PKZB is that it contains an em- 
pirical parameter in its exchange term, which was fitted 
to molecular atomization energies. Consequently, PKZB 
produces too-long bond lengths and some inaccurate 
properties of hydrogen-bonded complexe s 19 ' 25 ! 26 . These 
failures may be attributed to the unbalanced descrip- 
tion of PKZB exchange and correlation for slowly- varying 
densities and one- or two-electron densities, which are the 
paradigms in condensed matter physics and in quantum 
chemistry, respectively. 

Starting with the PBE GGA, Tao, Perdew, Staroverov, 
and Scuseria (TPSS) 23 have constructed a nonempirical 
meta-GGA. While the formula for TPSS looks a little 
more complicated than the PKZB one, the guiding the- 
ory is simple. A sound meta-GGA should be able to 
describe well the paradigm densities of condensed mat- 
ter physics and quantum chemistry. By imposing correct 
constraints, a meta-GGA can be made accurate for di- 
verse systems of interest. The TPSS construction builds 



many additional correct constraint s 11 ' 2 into a meta- 
GGA, while retaining those that the GGA has already 
respected. As a result, this meta-GGA is uniformly ac- 
curate for various properties of diverse systems 2 ^, sug- 
gesting the correctness of the TPSS philosophy. 

For high-density (exchange dominated) systems such 
as atoms, the TPSS is remarkably accurat o 27 ' 29 . In 
the low-density or strong-interaction limit, TPSS recov- 
ers the PKZB correlation, which is accurate for spin- 
unpolarized densities. The relative spin polarization is 
defined as 



where n = n-j + n^. Since TPSS correlation strives to 
represent a self-correlation correction to PBE, it should 
not change PBE correlation for a system with delocal- 
ized electrons, whether spin-polarized or not. This re- 
quirement is satisfied by TPSS through design for a spin- 
unpolarized jelliumii. However, we never impose this re- 
quirement for a spin-polarized density. Instead we first 
scale to the low-density limi t 30 ' 31 ' 32 ' 33 and there require^! 
the exchange-correlation energy to be correctly indepen- 
dent of spin for a model uniformly spin-polarized one- 
electron Gaussian density with constant relative spin po- 
larization in the range < \Q\ < 0.7, like LSD and PBE 
and unlike PKZB. 

In this work, we investigate the spin dependence of 
TPSS correlation and find that it is nearly the same 
as that of the PBE for spin-polarized jellium gener- 
ated by magnetic fields, implying that TPSS does not 
alter the PBE correlation energy for a spin-polarized 
system with delocalized electrons. Since TPSS success- 
fully improves on LSD for spin-unpolarized jellium and 
has the proper spin dependence, we estimate the surface 
exchange-correlation energy and work function with the 
TPSS meta-GGA functional for a spin-polarized jellium 
in magnetic fields. An application to the infinite bar- 
rier model (IBM)2i of metal surfaces is a related test but 
for rapidly- varying densities. The other tests considered 
here are jellium spheres (which sample the surface and 
curvature energy) and the linear response of bulk jellium. 



II. DENSITY FUNCTIONAL 
APPROXIMATIONS 

Density functionals may be ordered by the "Jacob's 
ladder" 35 ' 36 of approximations, according to the type of 
their local ingredients, whether constructed nonempir- 
ically or empirically. Here we only focus on the all- 
purpose nonempirical functionals, LSD, PBE and TPSS, 
and the nearly nonempirical PKZB, but not the ones 
recently developed for solids and solid surfaces such as 
AM052£, Wu-CoheiA and PBEsoI 3 ^. The first three 
rungs of the ladder to be considered here are LSD, PBE, 
and TPSS. (Note that the GGA described in Ref. Il 
is constructed in part from TPSS by the approximation 
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qb ~ 2p/3 for slowly-varying densities, which is a misin- 
terpretation of a statement in Ref. 23; it is only for a = 1, 
as in the uniform gas, and not for a ~ 1, that q~b « 2p/3.) 

Because the exchange component of a density func- 
tional satisfies the spin scaling relation 4 ^ 



£ x [n T ,n;] = £ x [2n T ]/2 + £ x [2n ; ]/2, 



(3) 



where -E x [n] = £" x [n/2,n/2] and n(r) = nf(r) + n|(r), 
we only need an exchange functional -E x [n] of a spin- 
unpolarized system. An exchange functional also satisfies 
the uniform coordinate scaling requirement 41 - 



(4) 



where 



7 3 rt(7r) is the scaled density of n(r). These 



two constraints are the basic requirements of an exchange 
functional. 

For a spin-unpolarized (closed shell) system, the ex- 
change functionals of the first three rungs may be ex- 
pressed in the form 



E x [n] = / d 



f (n)F x , 



(5) 



where e x mf (n) 



-^(S^n) 1 / 3 is the exchange energy 
per particle of a uniform electron gas and F x is the ex- 
change enhancement factor showing the nonlocality^ 7 - 



1 + k — k/(1 + x/ k), 



(G) 



with k — 0.804 and x > 0. The order of the ladder rungs 
depends upon the choice of x in Eq. ©. We have x = 
for LSD, and x — [ip for PBE, where /i = 0.21951, and 



p = |Vn| 2 /[4(37r 



2^2/3^8/31 = s 2 



(7) 



is the square of the reduced density gradient s. For the 
TPSS and PKZB meta-GGAs, x is a function of the two 
variables p and z, where 



/t < 1 



(8) 



with r = 'Y^ la T a and with t w = i|Vn| 2 /n being the von 
Weizsacker kinetic energy density. In the uniform-gas 
limit, all the density functionals above the first rung cor- 
rectly reduce to LSD. This uniform-gas limit is the most 
important requirement^ 2 - for bulk solids and surfaces. 

Since the correlation component of a density functional 
does not have such a simple spin scaling relation as the 
exchange component, we have to build the spin depen- 
dence into the correlation part. The LSD correlation en- 
ergy has the form 



E 



LSDr 



[n hni ] = I d A r ne^nf.nj), 



(9) 



where e" nlf is the correlation energy per electron^ 3 - for a 
uniform electron gas. The PBE correlatio n 18 ! 44 is based 
on the LSD correlation 



E. 



PBE r 



nr.ru] = / d r nec ("T.ni.Vnr.VnJ, (10) 



where e^ BE correctly reduces to e" nlf in the uniform-gas 
limit. The TPSS correlation 23,2 - is constructed from the 
PBE correlation, 

£ c TPSS [n T ,nJ = [d 3 r ne^tnt.ni.Vnt.Vni.Tt.Ti), 



where 



(11) 



'[l + deT^^i^ /t)% (12) 



g TPSS _ £ rcvPKZBn i ^rcvPKZB ( JW /^3i 



The quantity e T c cvPKZB of Eq. O is the revised PKZB 
correlation 

rovPKZB PBE/ Y7 V7 \h i rU f CM W I \2i 

e c = e c (nr,nx, Vnr,Vnx)[l + C(C,f)(r /r) J 



[i+c(c,mr w /rrj2- 



(13) 



with eJT being 



max[e. 



PBE/ 



,0,Vn CT ,0), 



PBE 



(nr,ni, Vnr,Vnx)], 



(14) 



where C is the relative spin polarization defined as £(r) = 
[n T (r) - rj;(r))/rc(r)] and £ = | VC|/[2(3^ 2 n) 1 / 3 }S 3 . 1 2 I^. 
Here C(C,£) at ( = and £ = in Eq. ([T3]) is chosen 
so that, in the low-density or strong- interaction limit, 
TPSS correlation recovers PKZB correlation, which is 
accurat e - 27 ! 33 for spin-unpolarized densities. The param- 
eter d in Eq. ([T!?]) is chosen such that the self-interaction 
correction should not alter the surface PBE correlation 
energy for spin-unpolarized jellium with delocalized elec- 
trons. The natural construction of the spin-dependent 
C(CjC) would be to make the TPSS correlation remain 
the same as the PBE correlation for spin-polarized jel- 
lium in the achievable (and thus energetically impor- 
tant) range of the uniform bulk relative spin polarization 
< C ^ 0.7. However, C(£,£) is designed instead to 
make TPSS independent of spin in the low-density limit 
for the one-electron Gaussian density and other densi- 
ties with uniform £ in the range of < C 0.7 without 
changing other properties. We will show that these two 
different procedures are essentially equivalent in the con- 
struction of the spin-dependent parameter C(£,£). 



III. SPIN-POLARIZED JELLIUM AND 
KOHN-SHAM APPROACH 

In a semi-infinite jellium metal filling the half space 
x < 0, the uniform positive background density corre- 
sponding to the ion lattice may be written as 1 - 



n + (r) = n8(—x), 



(15) 



where 9(—x) is a step function and n is the bulk electron 
density. In the absence of an external electric field, the 
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TABLE I: The strength [ibBq of two model external magnetic fields of Eqs. (|31f) and (|32[) . work function W, and surface 
correlation and exchange-correlation energies of the planar jellium surface in LSD, PBE, and PKZB and TPSS, as functions of 
uniform bulk relative spin polarization £ at normal bulk densities r s = 2,4,6. ^bBo and W are in eV; surface energies are in 
erg/cm 2 . LSD orbitals and densities are used. (1 hartree = 27.21 eV; 1 hartree/bohr 2 = 1.557 x 10 6 erg/cm 2 ) 
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correlation potential, and i> ex t(r) is the external scalar 
potential 

« cxt (r) = -/dV^. (19) 

(7 = 4-1 corresponds to electrons with spin a parallel to 
B and —1 to electrons with spin a antiparallel to B. The 
electron density n a may be evaluated from the occupied 
Kohn-Sham orbitals via 

n a (r) = £hAf| 2 ^-<), (20) 

where eg = k£ 2 /2 with k% = {%i: 2 n a ) 1 / :i is the Fermi 
energy per electron of spin a. Since of Eq. ([TT|) 

depends upon the density n a of Eq. ([20]), Eqs. ((171)- PD]) 
must be solved self-consistently. 

Suppose the external magnetic field is uniform for 
x < (within the bulk metal). Then the bulk density is 



electron density n(r) is related to the background density 
via the charge neutralization condition^ 

J d 3 r[n(r) - n+(r)] = 0. (16) 

In the presence of an external magnetic field B(r), the 
self-consistent single-particle Kohn-Sham equation may 
be written as^ 

j-lv 2 + vl ff {v)^\v) = e^f(r), (17) 

where v"* * is the local effective potential given by 

<//(r) = u«t(r) + J d 3 r' - a ■ M B B(r) 

+ <c(r). (18) 

Here er • B = <jB with a = ±1, /i^ = eh/ 2m is the 
Bohr magneton, v° (r) = 5E XC / '8n a (r) is the exchange- 
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spin-polarized uniformly with position-independent rela- 
tive spin polarization. At the Fermi levels, the spin-up 
and spin-down electrons have the same chemical poten- 
tial 

MT = AHi ( 21 ) 

where 

\i a = SE/6n a (r) 

ST 

= t — f-r +Vext(r) +Ujj(r) -afj, B B+ < c (r). 
dn^r) 

(22) 

Here T s [n-|-,njJ is the Kohn-Sham kinetic energy evalu- 
ated via Eq. (|TJ) and uh is the Hartree potential which is 
given by the second term of Eq. (| 18[) - The bulk uniform 
relative spin polarization is determined by the strength 
fj-sB of the applied external magnetic field via the rela- 
tion 

fi B B = (e F - 4)/2+ [vl(r s ,0 - ?4(r s , C)]A (23) 

where r s — (3/47rn) 1 / 3 is the bulk Seitz radius. For con- 
venience, we put the spin-up and spin-down Fermi levels 
at zero energy, so that we have 

<»<" — scfr (24 » 

Because the positive background is uniform, the electron 
density is also uniform in the bulk. Kohn-Sham wave- 
functions as solutions of one-electron Kohn-Sham equa- 
tion (fT?)) are plane waves and v^^(-oo) = — e F . The 
electrostatic potential is thus given as 

w os (-oo) = -e F - v° c (n, C) + ct^bB, (25) 

where v cs = v ext + u H - 

Near or at the surface, we may writei^ 

tf(T)=M(x)e« h »>+ k '>\ (26) 

where k, k v , and k z are the magnitudes of the wave vec- 
tors along x, y, and z directions, respectively. The three- 
dimensional Kohn-Sham equation (|17[) reduces then to 
the one-dimensional one^ for ip^(x), 

(27) 

where e k = k 2 /2. Solving Eq. (JSJJ) for ip%(x) yields the 
density n CT via 

n a (x) = 3n a f dfc (1 - P)|V-fc(x)| 2 , (28) 
Jo 

where < k = k/k F < 1. In the vacuum, the density 
decays exponentially^, as in an atomic, 

^(x) -> e" *, (a; -»«>), (29) 



TABLE II: Spin dependences of the surface exchange- 
correlation energies (in units of erg/cm 2 ) of the planar jellium 
surface in PBE and TPSS with various bulk valence-densities. 
LSD orbitals and densities are used. 
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where a is a constant for a given n. 

Within LSD, the exchange-correlation potential may 
be evaluated as 

d 

< c = ^-[« e xc(«T,«i)]- (30) 

Following Monnier and Perdew^ for the treatment of the 
electrostatic potential and the self-consistency procedure 
(outlined in Appendix A of Ref. [8j), we solved the one- 
dimensional Kohn-Sham equation (|2"T|) sclf-consistently 
within LSD. 

In the present work, two external magnetic fields cou- 
pled to the electron spins are considered: 

(1) Uniform inside and zero outside the jellium edge, 

B(x) = B 9(-x); (31) 

(2) Uniform everywhere, 

B(x) = B . (32) 

Eq. (f3"Tj) was proposed by Kautz and Schwartz^ to sim- 
ulate, within the jellium model, the "internal" magnetic 
field near the surface of a ferromagnetic metal, while 
Eq. (|32p can be realized experimentally over a range of 
Bq. The magnitude of the external magnetic fields Bq 
may be found from Eq. (|23[) for a given 

The work function W^&& is an interesting quantity 
which can be measured experimentally. It is defined as 
the energy required to remove an electron from a bulk 
solid into the vacuum. Within the framework of Kohn- 
Sham density functional theory the work function is the 
difference of the Kohn-Sham single-particle energies of an 
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TABLE III: Spin dependences of the surface exchange- 
correlation energies of the planar jellium surface in PBE and 
TPSS, when the enhancement factor _F XC of Eq. ([36} is uni- 
formly scaled to the high-density (r s — ► 0) or exchange-only 
limit from normal bulk valence-densities. LSD orbitals and 
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1265 


1709 


4 141 


PBE 


12 


39 


73 


107 




TPSS 


13 


40 


74 


108 


6 15 


PBE 


2 


7 


13 


19 




TPSS 


3 


7 


14 


15 



electron at rest in the vacuum and an electron moving at 
the Fermi level in the bulk 8 , i.e., 

W = < s (oo)-[fcf/2 + < ff (-oo)]. (33) 

In the presence of an exernal magnetic field, the work 
function is given byit 

W = < ff (oo) - [kf/2 + < ff (-oo)] + a^ B B 
= [v cs (co) - v ea (-oc)] - v° c {n^,n{) - fcp 2 /2 
+(t^ B B, (34) 

where the second equality can be obtained by combining 
Eqs. (HHJ) and (gSJ). 



IV. SURFACE ENERGIES OF 
SPIN-POLARIZED JELLIUM 

A. Surface energies 

The surface energy of a solid is the energy required to 
split the solid per unit area of new surface formed. Here 
we only focus on the exchange-correlation component 
<7 XC of the total surface energy. The surface exchange- 
correlation energy is 

/oo 
dx n(x)[e xc (x) - e xc (-oo)], (35) 
-oo 

and may be decomposed as a sum of the exchange and 
correlation contributions <7 XC = er x + <j c . 



TABLE IV: Spin dependences of the surface exchange- 
correlation energies of the planar jellium surface in PBE and 
TPSS, when the enhancement factor F xc of Eq. ([36} is uni- 
formly scaled to the low-density (r s — > oo) limit from normal 
bulk valence-densities. LSD orbitals and densities are used. 
The (-dependence here arises not so much from F xc as from 
the ("-dependence of the surface density profile n(x), which 
for B(r) — Bo spreads out more as |£| increases, (erg/cm 2 ) 
















(r.,0)] 


r a 


(erg/cm 2 ) 




0.2 


0.4 


0.6 


0.8 








B(r)- 


= B 6{-x) 


1 




2 


6364 


PBE 


-17 


-77 


-215 


-504 








1 Q 


81 
01 


990 




4 


503 


PBE 


1 


3 


2 


-10 






TPSS 


1 


1 


-2 


-7 


6 


107 


PBE 


1 


3 


5 


3 






TPSS 


1 


2 


4 


6 








B( 


r) = Bo 






2 


6364 


PBE 


477 


1461 


2352 


2714 






TPSS 


470 


1441 


2352 


2753 


4 


503 


PBE 


26 


86 


149 


180 






TPSS 


25 


83 


148 


193 


6 


107 


PBE 


6 


18 


31 


37 






TPSS 


6 


17 


31 


42 



We evaluated the surface exchange and correlation 
energies of spin-polarized jellium produced by the two 
model external magnetic fields of Eqs. ([3"Tj) and (|32[) . 
The results are displayed in Table Q] The surface ex- 
change and correlation energies of spin-unpolarized jel- 
lium are also listed for comparison. In our calculations, 
we employed LSD orbitals and densities obtained by self- 
consistently solving the Kohn-Sham equation (|27|) using 
the LSD exchange-correlation potential. For a justifica- 
tion of this approach, see Ref. H3- The work functions 
shown in Table I are LSD values. 

Previous studies^ show that, like PKZB, TPSS suc- 
cessfully improves upon LSD in the surface exchange- 
correlation energy of a spin-unpolarized jellium, while 
PBE gives a correction of the wrong sign to LSD and 
thus underestimates this quantity. Figures 1 and 2 clearly 
show that, while the surface exchange energy of LSD is 
much more overestimated than those of the PBE GGA 
and TPSS meta-GGA, the surface exchange-correlation 
energies of these density functional are not very differ- 
ent from each other. Fig. 2 shows that the improvement 
of TPSS over LSD and PBE may be attributed to its 
recovery of the known correct gradient expansions of the 
exchang o 23 i 24 -correlatio n 48 i 49 i 50 energy. Pictured in Fig. 
2 is the dependence of the exchange-correlation enhance- 
ment factor, 

F xc = £x T c PSS (n t ,r li ,V flt ,V Jll ,r t ,r i )/ £ r f (n) I (36) 
upon position x. 
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FIG. 1: Exchange enhancement factor F x of Eq. ([5]) for LSD, 
PBE, and TPSS as functions of position x in units of 2-n /fcp 
relative to the jellium edge (x = 0) for spin-unpolarized jel- 
lium with r 3 — 2. 
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FIG. 2: Exchange-correlation enhancement factors F xc of 
Eq. (J35J) for LSD, PBE, and TPSS as functions of position 
x in units of 2n/kF relative to the jellium edge (x = 0) for 
spin-unpolarized jellium with r 3 = 2. The same (LSD) surface 
density profile n(x) is assumed for all three functionals. Thus 
the differences in F xc (x) shown here determine the differences 
in surface exchange-correlation energy <r xc seen in the £ = 
row of Table I. The higher is F xc at a given x, the lower is 
cr xc . But the large- a; tail region of low density is much less 
important than the region close to the edge of the positive 
background at x = 0. 



Table I suggests that the work function and surface 
exchange-correlation energy of a metal could be slightly 
reduced by an increase in bulk spin polarization due to 
ferromagnetism [B = Bq6(— x)], but could be strongly 
increased by an increase in spin polarization due to an 
external uniform magnetic field. For the work function 
of the model ferromagnet, our results agree qualitatively 
with those of Kautz and ScwartzS (who however used 
random phase approximation (RPA) input to their LSD 



FIG. 3: Comparison of TPSS exchange-correlation enhance- 
ment factor _F XC of Eq. (|36|l as a function of x in units of 2n /k-p 
relative to the jellium edge (x = 0), for spin-unpolarized jel- 
lium at r a = 4, and for spin-polarized jellium at r s — 4 and 
C = 0.4 produced by two external magnetic fields of Eqs. (PTj) 
and (EU). 
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FIG. 4: Relative spin polarization £ = (n-f — ny)/n as a func- 
tion of x in units of 2n /fcp relative to the jellium edge (x — 0) 
for spin-polarized jellium at r s — 4 and £ = 0.4 produced by 
two external magnetic fields of Eqs. (13 1 p and (|32[) . 



calculation) . 

The local spin polarization varies from its bulk 
value at x = — oo to a limiting tail value C(°°) at x = 
+oo. For B = B 6(-x), ((oo) = 0, but for B = B , 
C(oo) = 1. These different limits are reflected in Fig. 
3, a plot of exchange-correlation enhancement factors vs. 
position x, and are presented directly in Fig. 4. 



B. Spin dependences of PBE GGA and TPSS 
meta-GGA 

To examine the spin dependences of the PBE GGA and 
TPSS meta-GGA for a spin-polarized jellium, we define 
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TABLE V: Surface exchange-correlation energy of jellium 
in the infinite barrier model. March's unit 3e 2 n/47r = 
(88738/r?) (erg/cm 2 ) is used^. The RPA+ value for r s = 
is from Ref. [5a . 



r s 


RPA + 


LSD 


PBE 


PKZB 


TPSS 





0.0714 


0.1108 


0.0452 


0.0501 


0.0516 


2 


0.1289 


0.1222 


0.1035 


0.1079 


0.1098 


4 


0.1364 


0.1296 


0.1190 


0.1235 


0.1255 


6 


0.1393 


0.1350 


0.1286 


0.1335 


0.1354 


oo 




0.2157 


0.2313 


0.2588 


0.2603 



the surface exchange-correlation spin-polarization energy 
as 

Acr xc = <7 xc (r s , C) - cr xc (r s , 0). (37) 

This quantity Act xc can tell us how the surface exchange- 
correlation energy of a spin-polarized jellium changes 
when ( changes. Tables UH IIII1 and IIVI show the com- 
parison of the spin dependences of PBE and TPSS for 
normal bulk valence densities and when the enhancement 
factor F xc of Eq. (f3"tj|) is uniformly scaled to the high- 
density (r s — > or exchange-only) and low-densit y 27 i 33 
(r s — > oo or strong-interaction) limits, respectively. From 
Tables HI Ell and [IV] we see that TPSS has nearly the 
same dependence as PBE in every case we examined here, 
suggesting that our earlier procedure to construct the 
spin-dependence of C(C,£) m Eq. (flU|) is right. 

The PBE GGA and TPSS meta-GGA describe mag- 
netism very similarly at jellium surfaces. Whether this 
will remain true for real systems remains to be deter- 
mined. We only know of a study of the ground-state spins 
of iron complexes^, which compares PBE with TPSSh 19 
(a hybrid of TPSS with 10% exact exchange), and an- 
other study of iron complexes^ which seems to show 
similar energy gaps between high- and low-spin states 
from PBE and TPSS. A private communication^ from 
one of the authors of Ref. [5l| shows that PBE and TPSS 
give similar energy differences among three spin states in 
each of seven iron complexes. 



V. INFINITE BARRIER MODEL OF THE 
JELLIUM SURFACE 

The earliest surface model of a metal is the infinite 
barrier model proposed by Bardeen^, in which the sur- 
face density profile is that of noninteracting free elec- 
trons in the presence of a hard wall. This is an oversim- 
plified surface model of a jellium metal, because many 
properties of this model surface can not be transferred 
to real metal surfaces. However, it may serve as an 
ideal model system with rapidly varying densities and 
can be employed to test density functionals, because the 
exact electron density, the kinetic energy density, and the 
conventional exchange energy per electron of this model 
system are analytically know n 34 i 55 . Here, the surface 




-0.074-1 , 1 , 1 , 1 , 1 , 1 , 

20 40 60 80 100 120 

number of electrons 

FIG. 5: Total energies per electron of jellium spheres for r s = 
4. The effect of our fixed-node correction is visible for magic 
clusters N = 2, 8, 18, 20, 34, 40, 58, 92, and 106. PBE is 
too low, but TPSS and especially LSD are good for TV > 2. 
LSD is best for N > 2, but worst for N = 2. Although TPSS 
surface energies are slightly higher than LSD values, TPSS 
curvature energies (Fig. 6 and Table IX) are lower, leading 
to a slightly lower TPSS total energy for these spheres. All 
curves tend as N — > oo to —0.0774 hartree. 



exchange-correlation energies of the LSD, PBE, PKZB 
and TPSS are evaluated for normal bulk valence densi- 
ties, and when the enhancement factor F xc of Eq. (|36|) 
is uniformly scaled to the high-density or exchange-only 
and low-density or strong- interaction limits. The results 
are shown in Table fVl 

RPA+ is a sophisticated approximation involving the 
full random phase approximation (RPA) plus a GGA for 
the short-range correction to RPA. The RPA+ values in 
Table [V] are evaluated from^ 

cr^ A+ = 0.0714(1 + 3.451r a )/(l + 1.688r s ), (38) 

which is a fit to the RPA+ values of Yan, Perdew, and 
Kurth£2 at r s = 0, 2.07, 4 and 6. The RPA+ values 
are exact in the exchange-only or r s — > limit, and we 
take them to be nearly exact for all rs. We see from 
Table [V] that the order of accuracy of these functionals 
for normal bulk valence densities is er P c BE < °xc KZB ~ 
°xc PSS < Cm ! while, in the high-density or exchange- 
only limit, we have ct lsd < cr PBE < ct pkzb w aj pss . 



VI. ENERGIES OF JELLIUM SPHERES 

In the spherical jellium model a positive charge back- 
ground is contained inside a sphere of radius 

R = r s N 1/3 . 
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FIG. 6: Deviation from LSD of PBE, PKZB and TPSS en- 
ergies for jellium spheres with r s — 3.93. The full lines are 
parabolas fitted to the LDM via Eq. (|42[) , as explained in the 
text. The open circles are input values for the magic clusters 
respecting the "Aufbau" principle. The second derivative of 
this curve at 7V~ 1//3 = is proportional to the difference of 
curvature energies between the given functional and LSD. 



The potential due to the positive background charge is 

(r < R) 



v + (r) 



N 
' 2R 



(39) 



Calculation of the exchange and correlation energies was 
done in an a posteriori process using LSD densities as 
input, as explained in Ref. [U By solving the many- 
electron problem using the Kohn-Sham approach we may 
obtain a series of single particle levels. For metallic den- 
sities the energy ordering of these states is Is, lp, Id, 2s, 
If, 2p, lg, 2d, lh, 3s, 2f, ... The filling of shells yields 
special stability for the so-called magic clusters. 

As in Ref. |47|, we calculated the energies of jellium 
spheres for some magic numbers (N = 2, 8, 18, 20, 34, 
40, 58, 92 and 106) with various densities r s — 2.07, 3.25, 
4.00 and 5.62. We chose these magic numbers and densi- 
ties due to the availability of the respective results in Dif- 
fusion Quantum Monte Carlo calculations by Sottile and 
Ballone- 9 , with which we want to compare. These Dif- 
fusion Quantum Monte Carlo (DMC) results on jellium 
spheres are supposed to be the most accurate available 
for the systems under study. However, Sottile-Ballone 
values— are affected by a systematic error in the corre- 
lation energies due to their fixed-node assumption. 

For N — > oo we get the limit of the uniform electron 
gas for which the fixed-node error may be estimated as 
the difference between the fixed-node calculation done 
by Ortiz-Ballone^ and the released-node calculation of 
Ceperley-Alder^i. This error estimate for the uniform 
electron gas has been presented in Table VII of Ref. [47]. 
For a sphere with N — 2 electrons, there is no node 
in the space factor of the ground-state wavefunction so 
the fixed- node error is absent. Using the limits N = oo 
and N — 2 we may interpolate the fixed- node error of the 
intermediate spheres by multiplying the correction to the 
correlation energy in the uniform electron gas by a factor 
suggested by the Liquid Drop Model§2, (Eq. (|4T]) below): 



Ae c (iV) = Ae™ if 



N 



1/3' 



(40) 



(r > R). 



The uniform electron gas correction Ae"™ 1 ^ is equal to 
the Ae° s of Table VII of Ref. E3- 

A better comparison of the energetics produced by the 
various density functionals can now be done. The ef- 
fect of correction is shown in Fig. [5] for a single density 
(r s =4.0 bohr). The error of a density functional is the 
difference between the exchange-correlation energy and 
the corresponding corrected DMC value. Table |VT1 shows 
the errors of the different density functionals for the total 
energy per electron, and for the five indicated densities, 
averaged over the magic closed-shell clusters in the range 
2 < N < 106. 

Table rvLTl displays the relative errors in the correlation 
energies, again averaged over the closed-shell spheres. 
The improvement of all functionals with respect to LSD 
is clear and the TPSS shows the smallest deviation in 
correlation energy. 

Surface and curvature energies are relevant^ not only 
to clusters and voids, but even to cohesive energies 
and monovacancy formation energies. Adopting the 
same fitting procedure for extracting Liquid Drop Model 
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TABLE VI: Mean absolute deviations from corrected fixed- 
node DMC values [Ref. [5t| and Eq. (|4Up ] of the total energies 
per electron of jellium spheres in various density functional 
approaches. The values are averages over the magic clusters 
N = 2, 8, 18, 20, 34, 40, 58, 92, and 106. For individual N at 
r s — 4, see Fig. 5. 

\{E-E UMV )/N\ (hartree) 



r a 


LSD 


PBE 


PKZB 


TPSS 


1.00 


0.0040 


0.0015 


0.0015 


0.0007 


2.00 


0.0018 


0.0007 


0.0004 


0.0004 


3.25 


0.0007 


0.0005 


0.0002 


0.0004 


4.00 


0.0004 


0.0005 


0.0004 


0.0005 


5.62 


0.0005 


0.0006 


0.0006 


0.0007 


average 


0.0015 


0.0008 


0.0006 


0.0006 



TABLE VII: Average relative deviations of the correlation 
energy of jellium spheres, in various density functional ap- 
proaches, from corrected DMC values [Ref. [59| and Eq. (|40p ] . 
Averages were taken over magic clusters: N = 2, 8, 18, 20, 
34, 40, 58, 92, and 106. 







(E c -E» MO )/E» MU 




r s 


LSD 


PBE 


PKZB 


TPSS 


1.00 


40.3% 


6.7% 


7.4% 


5.9% 


2.00 


34.0% 


7.7% 


7.8% 


6.5% 


3.25 


29.7% 


7.3% 


6.9% 


5.8% 


4.00 


27.2% 


6.4% 


5.9% 


4.9% 


5.62 


26.4% 


7.3% 


6.4% 


5.5% 


average 


31.5% 


7.1% 


6.9% 


5.7% 



(LDM)£2 parameters from jellium spheres as in Ref. H3, 
we use the following equation for the energy of a neutral 
jellium cluster with N valence electrons 

■pLDM 

— = e mif + AnrlaN- 1 ^ + 2nr sl N- 2 / 3 , (41) 

where o and 7 describe the surface and curvature ener- 
gies respectively, and e uni f = (iirr^/fya is the energy 
per electron of the uniform electron gas, with a being 
its energy per volume. This model neglects the quan- 
tum oscillations in the energy due to the shell structure. 
For the sequence of closed-shell clusters, the oscillation 
is presumably the same in LSD as at any higher level 
of theory, so the difference cancels out. (More gener- 
ally, [E — Els~d)/N in finite systems can be extrapolated 
smoothly to infinite size^.) Thus the LDM equation for 
the closed-shell clusters, including the smaller ones, is 
better written as the difference to LSD 47 : 

K _ E LSD _ ( unif _ unif x 

+47rr s 2 (a - + 2 ^( 7 _ 

As the functional PBE, PKZB and TPSS have the same 
parametrization (PW92)^ in the limit of the uniform 
electron gas, the first term in the right-hand side is zero. 

We used surface energies calculated by a planar surface 
code of Monnier and Perdew 8 and reported in Table fVTlTl 



TABLE VIII: Surface exchange-correlation energies of jel- 
lium calculated with a planar surface code. Energies are in 
erg/cm 2 . Only deviations from LSD are relevant to fitting 
Eq. UHJ. 



r a 


LSD 


PBE 


PKZB 


TPSS 


2.07 


2961 


2881 


3002 


2985 


2.65 


1204 


1167 


1221 


1215 


3.24 


575.1 


555.9 


583.1 


581.8 


3.93 


279.8 


269.9 


283.8 


284.0 


5.62 


69.89 


67.27 


71.15 


71.90 



TABLE IX: Curvature total energies of jellium (in units of 
millihartree/bohr) of jellium. Only deviations from LSD are 
relevant to fitting Eq. (l42jl . 



r s 


LSD 


PBE 


PKZB 


TPSS 


2.07 


1.830 


1.494 


0.988 


1.063 


2.65 


1.044 


0.885 


0.548 


0.609 


3.24 


0.635 


0.546 


0.312 


0.358 


3.93 


0.369 


0.318 


0.156 


0.189 


5.62 


0.180 


0.161 


0.082 


0.097 



Thus only the curvature-energy term in Eq. (|42[) needs to 
be fitted. To perform this fit we extended the calculation 
of jellium spheres up to N = 1100 (except in the case 
of r s =5.62, where we only reached up to N =748). We 
checked the "Aufbau" principle, i.e., the highest occupied 
Kohn-Sham orbital should have an energy lower than the 
lowest unnoccupied orbital, and we only took the clusters 
obeying that principle. An example of such fitting is 
plotted in Fig. [Sj 

Using the resulting fit to curvature-energy differences 
and the LSD curvature energies given by Ziesche et al£^, 
we show in Tabic IIXI our curvature energies for several 
densities. The energies are in fact very similar to a pre- 
vious calculation 4 !, which used only a smaller number of 
spheres and did not restrict the surface energy term. The 
curvature energies of TPSS are close to those of PKZB, 
as expected. The PKZB curvature energies are smaller 
than those of LSD and PBE, as predicted in Ref. H3. 



VII. LINEAR DENSITY RESPONSE AND 
CHARGE DENSITY WAVES 

There are several reasons why the linear response of 
the homogenous electron gas to an external perturbation 
is of interest. First, it is an important step towards a 
qualitative understanding of the electronic structure of 
the simple metals, via pseudopotential perturbation the- 
ory^. Second, comparing the response obtained from ap- 
proximate density functionals to the one obtained from 
Quantum Monte Carlo calculations can serve as a test of 
different density functionals. Third, the linear response 
relations can be used to study the instability of the uni- 
form phase of jellium against the formation of a charge- 
density wave§£. 
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FIG. 7: The exchange-only response function 7 X for bulk jel- 
lium in the TPSS (full line), PBE (short dashed) and LSD 
(dotted) functionals. The long dashed line shows the exact- 
exchange only results from Ref. l69L See text for discussion. 




0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 

y=k/(2 k F ) 

FIG. 8: The exchange-correlation response function 7 XC for 
bulk jellium as obtained from the TPSS functional (full line), 
the LSD and PBE functionals (long dashed), the Richardson- 
Ashcroft (RA) approximation^- (short dashed) and the Quan- 
tum Monte Carlo calculations of Ref. [z3 (crosses) for r s = 2. 



In the present work, we calculate the linear response 
for LSD 43 , PBE GGAi£, and the PKZB^ and TPSS2 3 . 
meta-GGAs. The quantity of interest is the response 
function 7xc(q) defined by the linear response relation 



Sv xc (k) = --rjTxc ( — ) Sn(k). 

rbT? 



k 
2k F 



(43) 



This equation relates the Fourier component 5n(k) of a 
density perturbation to the Fourier component Sv xc (k) of 
the exchange-correlation potential that results from the 
perturbed density. We calculated 7 XC separately for the 
exchange and correlation parts of the different functionals 
by inserting the perturbed density n(r) = no+n^ cos(kr), 
where n = k F /(3ir 2 ) and % < iioj m to the energy func- 
tional. The resulting expression can be rewritten as a 



FIG. 9: Same as Fig. 8 but for a density with r s = 5. 



power series in n^/no and 7 is obtained by multiplying 

the second-order coefficient of this expansion by — ^— s -. 
(The factor 2 takes care of the 1/2 in the Taylor expan- 
sion.) This procedure straightforwardly yields 7x SD = 1 



and 7x BE = 1 + f^y 2 , where y = jj^ and fi = /37r 2 / 3 , 
with (3 — 0.066725 being the coefficient of the second- 
order gradient term in the gradient expansion of the cor- 
relation energy in the high-density limit^. 

For the meta-GGA functionals the situation is more 
complicated, because they also depend on the kinetic en- 
ergy density r, which must be calculated from the or- 
bitals. We can, however, obtain the linear response with- 
out knowing the orbitals in certain limiting cases by using 
the appropriate gradient expansions of the kinetic energy 
density. For a slowly- varying perturbation, i.e., k — > 0, 
the gradient expansion reads^i 



1 



V 4 n 



360 (37r 2 ) 2 / 3 n 2 / 3 ' 



(44) 



We here retained only terms up to order nk since higher 
orders will not contribute to the linear response. In- 
serting this ex pre ssion into the TPSS energy functional, 
Eq.(3) of Ref. I23L and again expanding into a series in 
n k/%j we obtain after some algebra the slowly- varying 
limit of the TPSS linear response function 



Tx.slow 1 



73 

225 



V 



146 
3375 



y 6 + 0(y s ). (45) 



In the rapidly- varying limit, i.e., k — * 00 and <C 1, 
the kinetic energy density can be expanded as^ 



T v 3.(3^2)2/3 6/3 + ILY^! + cV V (46) 
10 v 1 8 n ' 



To the best of our knowledge, no value for the coefficient 
c has been given in the past. We argue that c = for 
the following reason: It is known that the von Weizsacker 
kinetic energy density is a rigorous lower bound for the 
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kinetic energy density, 



r > 



1 |Vn| 2 



(47) 



But for c ^ and any non-vanishing amplitude of the 
perturbation, one can choose a wavevector to make the 
V 2 -term in Eq. (|4"6"| arbitrarily large and negative for cer- 
tain points in space. Thus, the only way to avoid the 
violation of the inequality (|4T[) is to simply require c = 0. 
Using Eq. (|46[) with c = yields the rapidly- varying limit 
of TPSS linear response, 



TPSS 
ix, rapid 



i 



9 



y 



(48) 



Finally interpolating between the two limiting cases of 
Eq. pS)) and Eq. (|4"8")) with the Pade approximant 



7 p,tpss = j 



9 V 



73 4 
225^ 



(1 



(49) 



we obtain an expression that can be expected to be very 
close to the exact TPSS response 7j PSS for all k. Going 
through the same procedure for the PKZB meta-GGA 
we confirmed that PKZB and TPSS have the same linear 
response. 

In Fig. 7 we plotted the exchange-only response func- 
tions for LSD, PBE GGA, and TPSS meta-GGA, to- 
gether with the exact exchange-only response from Ref. 
l69l The TPSS exchange-only response is extremely close 
to the exact exchange-only response up to y sa 0.6, and 
both PBE and TPSS provide reasonable approximations 
up to y w 1. Only for wavevectors with a magnitude of 
more than twice the Fermi wavenumber do differences be- 
come pronounced since the semi-local functionals do not 
recover the abrupt drop of the exact exchange response 
beyond y s» 1. 

To calculate the correlation contribution to the LSD 
linear response, we proceed slightly differently and di- 
rectly evaluate 7c LSD = ~{k 2 /i:)5 2 E^ /5n(r)5n(r') 
with the Perdew-Wang expression^ 3 -. Obviously, 
S 2 E^ SD /5n{r)Sn(r') = [2de^ ni /dn + nd 2 e^ ni /dn 2 ] S(r - 
r'), and some algebra yields 



4nri 



On 



-A 



Cl<72 



r s c 3 qi 



2ai In 1 



(50) 



and 



2 ,uni 



d 2 e 



dn 2 



1G 

81 



7T r s [aq 2 + Aaiqiq 2 c 3 



2ciq 2 c 3 



+ cic 3 (-/3i + 3p 3 r s + 8/3 4 ^ /2 )/V? 
/[8^(l + l/c 2 ) 2 ] 
Aciq 2 



Here qi = Pi 

2/3 2 Jr; + 3(3 3 r 



l^/r s qic 3 



+ P2 



2Aa x ln(l + 1/ca) 



(51) 



A a 3 / 2 

4/3 4 r s , ci 



a 3/2 

p 3 ?v 

= 1 + a\r s , c 2 = 2aq\, and 



far 2 , q 2 = Pi + 



c 3 = 1 + cgjwith A, ai, and Pi— Pi being the parameters 
from Ref. |43| . 

The contribution of the correlation part to the response 
for the PBE GGA and the two MGGAs is obtained anal- 
ogously to the calculations for the exchange functionals. 
In all three cases we obtain the same result, 

7 PBE = 7 PKZB = 7 TPSS = 7 LSD _ 3^ (g2) 

where p is the same parameter as in the PBE exchange 
response function (see above). 

In Figs. 8 and 9 we compare the response as obtained 
from the density functionals to the Quantum Monte 
Carlo results of Ref. [73 and the results of Ref. [7l| for 
two different densities. LSD and PBE provide a satisfy- 
ing average 7 for y < 1. The TPSS response shows more 
structure than LSD and PBE and is also in satisfactory 
agreement with the QMC results. Only for y > 1.1 does 
the TPSS response move outside of the QMC error bars. 

Finally, following Ref. [66J we calculated the instability 
of jellium against the formation of a charge density wave. 
The Fourier components of the self-consistent potential 
are related to the ones of the external potential by 



«(k) 



1 



2/ 2 7x 



Uext(k), (53) 



where x is the Lindhard response function. For vanish- 
ing « ext (k), the left-hand side of this equation can only be 
nonzero if the denominator on the right-hand side van- 
ishes. Thus, for each value of y we numerically search for 
the density (i.e., k-p) that sets the denominator to zero. 
The largest value of k-p (for all y) found in that way 
marks the onset of jellium instability. In this way we 
confirmed that, for exchange and correlation combined, 
TPSS does not alter much the prediction of LSD for the 
onset of jellium instability, which occurs for fcp < 0.06 
a.u. (r s > 30). 



VIII. CONCLUSIONS 

In conclusion, we have calculated and compared the 
jellium surface exchange-correlation energies of the PBE 
GGA and the TPSS meta-GGA and of these two func- 
tionals when uniformly scaled to the high- and low- 
density limits for the normal bulk valence densities in 
magnetic fields. In all the cases, the fairly good agree- 
ment of the PBE GGA with the TPSS meta-GGA shows 
that the TPSS meta-GGA indeed represents the self- 
correlation correction of the PBE GGA. We have further 
found that the "internal" magnetic field of Eq. (f3Tj) and 
the external uniform field of Eq. (|32|) are typically oppo- 
site to each other in their effects on the work function 
and surface exchange-correlation energy of jellium. 

We have also calculated the energies of jellium spheres 
with LSD, PBE GGA, and TPSS and its predecessor 
PKZB meta-GGAs. Typically, while PBE energies are 
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too low for spheres with more than about two electrons, 
LSD and TPSS are accurate there, up to 106 electrons. 
Curvature energies are reduced substantially as we pass 
from LSD to PBE to TPSS. Finally, we have shown that 
the linear response of bulk jellium (to perturbations with 
wavevectors less than twice the Fermi wavevector) is rea- 
sonably described by all the functionals considered here. 

As we climb the ladder of nonempirical density 
functional approximations from LSD to GGA to 
meta-GGA, there is a steady and dramatic improve- 
ment in atomization energie s 19 ! 72 . Surface energies 
worsei g£JJ a i2 1 i3 1 i4j£ 1 i6 1 ^ f rom LSD to PBE GGA (and 

other popular GGAs), due to an imperfect error cancel- 
lation between exchange and correlation, but this can be 
corrected in any of three ways: (1) by transferring 7 -^ the 
needed correction from jellium to real systems, (2) by us- 
ing GGAs designed specifically for solids (and not for free 
atoms )2L2£i22, or (3) by climbing up further to the TPSS 
meta-GGA. The third way adds little in computational 
cos t 19 ' 72 , even at full selfconsistency, and seems worthy 



of further testing and possible refinement. 

The jellium model itself remains useful as a testing 
ground for density functionals. Although some of its 
properties become unphysical as one moves away from 
the bulk density at which jellium is stable (r s « 4, 
roughly the valence density of sodium), this problem 
can also be fixed inexpensively via the stabilized jellium 
models. 
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